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Abstract 

We discuss a new method for obtaining the WKB approximation to the Dirac equation with a 
scalar potential and a time-like vector potential. We use the WKB solutions to investigate the 
scaling behavior of a confining model for quark-hadron duality. In this model, a light quark is 
bound to a heavy di-quark by a linear scalar potential. Absorption of virtual photons promotes 
the quark to bound states. The analog of the parton model for this case is for a virtual photon to 
eject the bound, ground-state quark directly into free continuum states. We compare the scaling 
limits of the response functions for these two transitions. 

PACS numbers: 12.40.Nn, 12.39.Ki, 13.60.Hb 
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I. INTRODUCTION 



Quark-hadron duality implies that in certain kinematic regimes, properly averaged 
hadronic observables can be described by a perturbative QCD (pQCD) calculation. This 
version of duality is highly relevant as perturbative QCD calculations can be performed. 
Using duality, these pQCD calculations can then be related to averaged data taken in the 
resonance region. Duality was first observed more than 30 years ago experimentally by 
Bloom and Oilman in inclusive inelastic electron scattering and has since then been 



observed in a variety of reactions: e hadrons is an example for dualit};^we 
the textbooks, the semileptonic decay of heavy mesons is another example 



inow from 
duality 



is considered in the analysis of heavy ion reactions j^, and forms the basis for using QCD 
sum rules IgI . Recently, duality has been observed to high precision and down to rather low 

pnnn 

momentum transfers in electron scattering at Jefferson Lab [7|, |8|, |9|, 110|| . Duality in spin ob- 
servables is curre ntly pr obed by several experiments, both at Jefferson Lab and at DESY in 



IS currently pr( 



Germany |lll.ll2l.ll3Lll4|. In addition to the "classical" examples and applications of duality, 
duality ideas are applied in new areas, too. For neutrino scattering, the beam energies are 
not well known, and an averaging will thus take place almost automatically. The application 
of duality is discussed for several planned neutrino experiments, see e.g. jl5|, and duality 
ideas have been applied in 16] to nucleon/nuclear duality in neutrino scattering. There is 



also interest in duality in parity violation experiments 



parton distributions 



cyin par 



17|, and with regard to generalized 



A very local version of duality - assuming that it holds for 



20 



2l| to extract information on structure functions 



just one resonance - has been used in 

at XBj 1 in the scaling limit from form factor data. These ideas were also applied to 



neutrino- nucleon scattering 2l|. Duality ideas might also be useful for pion photoproduc- 
tion ^22]. Duality is a major point in the 12 GeV upgrade of CEBAF at Jefferson Lab {2^. 
A recent review of quark-hadron duality can be found in [241] . 

Apart from being interesting all by itself, quark-hadron duality is an important tool for 
studying kinematic areas that cannot be accessed in the deep inelastic regime: measurements 
at large values of xbj in the resonance region typically have much higher count rates than 
measurements at the same xbj in the deep inelastic regime, which requires very high four- 
momentum transfers Q^. Application of a proper averaging procedure to the resonance data 
may allow the extraction of deep inelastic information, e.g. in the case of the polarization 
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asymmetry Ai of the neutron for xbj 1- Before duality averaging procedures can be 
applied safely, we require a thorough understanding of where duality holds and how exact 



it is. Therefore, duality has been studied intensively by theorists 
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37 



38, 



during 



39, 



the past couple of 
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45 



. Most of the theoretical studies focus on duality in electron scattering, due to 
the large experimental program. The model we will discuss here is for electron scattering, 
too. 

Many recent papers have been devoted to modeling quark-hadron duality in simple, fully 

■S, to g ain a better understanding of the conditions under which 



solvable relativistic mode 
duality works |3, 27, 33, 



34 



35, 



40 



4i 



43l |. The idea of modeling duality is to capture 



just the essential physical conditions of this rather complex phenomenon. Typically, these 
basic requirements for a model are imposed: one requires a relativistic description of confined 
valence quarks, and one treats the hadrons in the infinitely narrow resonance approximation. 

In these papers, one important point was raised and discussed that is interesting not just 
for duality, but in general: do the scaling curves obtained assuming outgoing plane waves 
agree with the scaling curves obtained when we assume final state interactions? 

The general approach is to model the perturbative QCD (pQCD) picture by considering 
a quark bound within a potential in the initial state, and after the interaction with a virtual 
photon, the quark is considered "free" and the potential set to zero. This "bound-free" 
transition corresponds to a plane wave impulse approximation (PWIA). This transition is 
compared with a hadronic picture: the initial state is the same - the quark is bound in a 
potential - but after the absorption of the photon, it is in an excited, but still bound, state. 
This "bound-bound" transition corresponds to the excitation of resonances. 

In order to reproduce duality, a model must fulfill several conditions that are observed 
experimentally: first, the bound-free transition, corresponding to pQCD, must scale for large 
momentum transfers. Second, the bound-bound transition must scale, and the bound-bound 
and bound-free scaling functions have to agree so that the third condition can be fulfilled: 
at low momentum transfers, the bound-bound results should oscillate around the bound-free 
scaling curve. Note that in all models currently proposed, there are no gluons included, and 
therefore neither radiative corrections nor evolution of scaling curves are present. 



In our recent papers 



40 



4l|, we could show analytically that the two different scaling 



curves, found for the bound-free and bound-bound transition, do coincide. In these papers. 
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we considered a model where all particles were treated as scalars j40|, and a model where 



0, 



on 



J the quarks were treated as scalars, while electrons and photons had their proper spin 
41|. The treatment of the quarks as scalars considerably simplified our calculations, as the 
resulting Klein-Gordon equation could be solved analytically, by recognizing that it could be 
rewritten to resemble a Schrodinger equation. Once we introduced the proper spin for quarks 



42, |43[, we solved the Dirac equation numerically, using a Runge-Kutta-Fehlberg (RKF) 
algorithm [s^. The numerical accuracy we achieved this way restricted us to momentum 
transfers below q < 12GeV, a value where the calculated response functions had not fully 
converged to their scaling value yet. Even though we solved only for eigenenergies of up to 
12GeV, we found roughly 24,000 states below that energy. As for any relativistic problem, 
the density of energy states rapidly increases with the energy. 

Thus, we could not determine if the bound-bound and bound-free scaling curves were 
going to coincide, as for the other, simplified models we investigated previously. 

In a recent series of papers by Paris et al. similar models were considered, and 

solved numerically, and the result found there was a pronounced discrepancy between the 
bound-free and bound-bound transition. This led the authors to question the interpretation 
of the scaling curves extracted from deep inelastic scattering (DIS). In DIS, it is generally 
assumed that final state interactions are negligible, and that the scaling curves can be 
interpreted in terms of parton distribution functions. 

In this paper, we investigate the scaling behavior of our model with Dirac quarks. We 
do this by employing the WKB method to solve the Dirac equation numerically for very 
high momentum transfers. This allows us to investigate the scaling behavior of the response 
functions in the bound-bound transition at the relevant high momentum transfers. The 
bound-free transition does not require any complicated final state calculations, and can be 
evaluated in a straightforward manner. 

This paper is organized as follows: first, we remind the reader of the general ideas un- 
derlying the WKB approximation jsil . and discuss the WKB for solving the Dirac equation. 
Then, we introduce our model and present numerical results obtained in the WKB approxi- 
mation. We compare them to the results obtained by explicitly solving the differential equa- 
tion at g = 10 GeV, the highest value accessible with the Runge-Kutta-Fehlberg method. 
After validating the WKB calculations, we proceed to investigate the scaling behavior of the 
bound-bound and bound-free transitions at large q, where convergence has set in. 
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II. THE WKB APPROXIMATION FOR THE DIRAC EQUATION 



The Dirac hamiltonian for a particle in a scalar filed S{r) and a time-like vector field 
V{r) is 

H = ci-p + (3{m + S{r)) + V{r). (1) 
The wave functions for the Dirac hamiltonian are 

where the yi^{Qr) are the usual spin spherical harmonics. Then angular momentum quan- 
tum numbers are / = —k, + 1 for k < and / = k, for k, > with I = —k for k < and 
I = K — 1 for K > 0. The reduced radial wave functions G{r) and F{r) are solutions to the 
coupled equations 

hG'{r) + h-G{r) = {m + S{r) - V{r) + E)F{r) (3) 
r 

hF'{r) -h-F{r) = {m + S{r) + V{r) - E)G{r) . (4) 
r 

The differential equation can be solved numerically, but above a certain energy, finding 



5G 



52l | becomes 



the numerical solution with a standard Runge-Kutta-Fehlberg method 
inaccurate. In the relativistic treatment of a potential, the energy eigenvalues of higher 
states lie closer and closer together. At some high enough energy, the separation between 
neighboring states becomes so small that some states may be missed with shooting methods. 
Also, one always has to calculate all states consecutively, instead of being able to calculate a 
state with certain, given quantum numbers. While an improved method for the integration 
of the differential equation may be applied, the situation lends itself to the application of a 
semiclassical approximation, the WKB method j^. 

In order to illustrate that the assumptions of the WKB approximation are perfectly 
reasonable, we show a plot of a one- dimensional harmonic oscillator potential in Fig. ^ 
together with the exact solution of the wave function for an energy E = 26huj. One can 
see that the wave function oscillates in the classically allowed region, where E > V, and 
strongly resembles a plane wave there. Outside of the classically allowed region, the wave 
function is damped and goes to zero. 

The WKB approximation can be applied when the potential has only a small variation 
over several wavelengths of the particle, as is the case in Fig. ^near the origin. This feature 
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FIG. 1: The potential of a harmonic oscillator is shown, together with the exact wave function for 
an energy oi E = 26fiLO. 

makes it quite useful for studying the highly energetic excited states of a particle bound 
in a potential. The WKB approximation assumes that in the classically allowed regions, 
the wave function can be approximated by a plane wave with a position dependent effective 
wave vector. As the potential is almost constant over a few wavelengths, the change of the 
wave vector is small compared with the value of the local wave vector. 

The WKB approximation can then be viewed as an expansion in derivatives of the effec- 
tive phase. In the case of a one- dimensional wave equation, this is equivalent to an expansion 
in h and is, therefore, often referred to as the semiclassical or quasiclassical approximation. 
In the case of the radial equation for a three-dimensional Schrodinger wave equation, care 
must be taken to correctly treat the centrifugal barrier term that arises from the angular 
momentum operator acting on the angular functions. This term is also important for the 
solution of the corresponding classical problem, and needs to be included in the leading- 
order effective potential for the WKB approximation, even though it will involve powers of 
h. The WKB approximation is then carried out by attaching an expansion parameter to the 
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derivative terms in the wave equation and then proceeding as in the one- dimensional case. 

To implement the WKB approximation for the Dirac equation the radial equations Q 
and (jH) are modified by making the substitution Tin ^ k and replacing the remaining factors 
of h by the expansion parameter r] to give 

r/G'(r) + ^G(r) = {m + S{r) -V{r) + E)F{r) (5) 
r]F'{r)-'^F{r) = {m + S{r) + V{r) - E)G{r) . (6) 

We now can parameterize the wave functions for positive energy solutions as 

G(r) = e^^^*") (7) 

and 

F(r) = B{r)ei^^''^ (8) 

where ^(r) is the local phase for a plane wave and B{r) allows for differences in the upper 
and lower component wave functions. Substituting these into (0) and ^ yields 

i^'(r) + ^ = {m + S{r)-V{r) + E)B{r) (9) 
r]B'{r) + iB{r)^'{r) - ^5(r) = (m + S{r) + V{r) - E) . (10) 

Equation Q can be solved to give 

B(,.) = . (11) 

Substituting (fTT|) into (fTUI) produces the second order differential equation for the phase 

ivCir) - ^'{rf - IT] ^'}''\ ~ ^'j-""} e'(r) - ^ - (m + S{r)y + {V{r) - Ef 



+ ^'^'^ - ""^'^ ] = (12) 
'r\r m + S{r)-V{r) + E) ^ ^ 

The WKB approximation is obtained by expanding the phase function in powers of rj as 

oo 

e(r) = $^r/"a(r). (13) 

n=0 

Substituting (fT!^ into (fT^ gives 

n=l n=0 m=Q n=l 

. _.o K (\ S'{r\-V'{r 



-- - (m + S{r)f + (\/(r) - Ef - vi- \- ^ ^c^^ t//\\ f ) = • (^4) 



Equating coefficients of like powers of t] gives 



<o\r) - ^ - (m + Sir)Y + (\/(r) - Ef = 



(15) 



for n = 0, and 



S'(r)-V'ir) nfl 



m + S{r) - V{r) + E 



S'{r) - V'{r) 



r \r m + S{r) — V{r) + E 



(16) 



for n = 1. 

Equation (|T5|l can be solved to give 



e^(r) = ±W (y (r) - Ey - (m + S{r)f 



(17) 



and (fT^ can be solved to give 

I Coir) i S'{r) - V'{r) 



where 



and 



S'{r) - V\r) 



2 Co{r) 2 m + S{r)-V{r) + E 2^^(r)r Vr m + S{r) - V{r) + E ^ 
t^{C,{r)) ± hir) , ' (18) 



/ci(r) 



1 K 



1 

- + 



S\r) - V\r) 



»(e^(r)) 



r m + S(r) — Vij) + E 
S\r) - V'{r) 



2kQ[r) r 

2 ^^(r) 2 m + S{r) - V{r) + E 



(19) 



d 



d 



H^/hi^) - — \n{^/m + Sir)-Vir) + E) . 



(20) 



Keeping terms to order r], and using (fTTj) and (fTHj) the phase function can be integrated 
to give 



e(r) 



rfr'(e;(r')+<(r')) 
dr' 

/T 
dr' {ko{r') + riki{r')) (21) 



±A;o(r') + irj-^ \n{Jko(r')) - ir]-^ \n{Jm + S(r') - V(r') + E) ± iqkAr') 
dr 



dr' 



The upper-component radial wave function can then be written as 



^. , ^.^m + S{r)-V{r) + E 
G(r) = rr^ 



±- I dr'{ko{r') + r]ki{r')) 

T] 



(22) 
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where is the normahzation constant. 

The lower-component wave function can be obtained directly by using (|TT|) to first order 
in 7] yielding 



F(r) = Af- 



1 



^/k^^m + S{r) - V{r) + E 
i 

.V 



K 



r 



X exp 



dr'{kQ{r') + r]ki{r')) 



(23) 



From this point we will set r/ = 1 for convenience. 

In constructing the above solution we have assumed that we are examining the solution in 
the classically allowed region where the energy E is greater than the effective potential. For 
the potentials which we will use, there is only one interval in r which is classically allowed. 
It is bounded by the classical turning points r± which are the solutions to 

fco(r±) = 0. (24) 

Quantization of the bound states is obtained by requiring that 

dr'{ko{r') + k,{r')) = (^n - vr , (25) 

where n > 1, as in the Schrodinger case. The solutions can be extended to the classically 
forbidden regions by analytic continuation. Since the k^lr) vanishes at the classical turning 
points, the local wavelength A(r) = 27r/A;o(r) becomes arbitrarily large near the classical 
turning points in contradiction to the basic assumption of the WKB approximation. The 
expansion in rj therefore does not converge in the vicinity of the classical turning points. 
Techniques for matching the wave functions at the turning points and replacing them with 
smooth approximate wave functions near the turning points are well described in most 
introductory graduate quantum mechanics texts. There is, however, one additional compli- 
cation in this case associated with the first-order phase ^i(r) which does not appear in the 
Schrodinger case. This requires us to modify the usual approach to approximating the wave 
functions near the classical turning points. The solution to this complication is discussed in 
Appendix A. 

We will need to construct the wave functions in three regions: the classically forbidden 
region where < r < r_ (Region I); the classically allowed region where r_ < r < r+ (Region 
II); and the classically forbidden region where r+ < r (Region III). The wave functions must 
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be matched at the boundaries. The procedure for the upper-component wave function is the 
same as for the Schrodinger case. Assume that we start by choosing the root in Eq. ()17|) 
where ^Q^r) = +kQ{r). The wave function is then obtained by choosing 



M=-e U 



and then defining the wave function as 



Gn{r) = G{r)+G*{r). 



This gives 



Gn{r) = N 



+ S{r) - V{r) + E 
7^ 



cos 



TT 



dr\k,{r') + h{r'))-- 



(26) 



(27) 



(28) 



Applying the same procedure to the lower component yields 

1 



Fair) 



X 



^fk^^m + S{r) - V{r) + E 



- %ar))) cos (^^^ dr'{ko{r') + k,{r')) - ^ 
=F [koir) + ki{r)) sin (^j^ dr'{ko{r') + ki{r')) - ^ 
In the classically forbidden regions the wave vectors become complex with 



(29) 



and 



-ikn(r) = kn(r) 



—ikiir) = ki(r) 



K 



— + {m + S{r)y - {V{r) - Ey 



K K S'(r) - V'(r) 



(30) 



(31) 



2A;o(r) r m + S{r) - V{r) + E _ 

The wave functions will be either growing exponentially or exponentially damped. Since we 
require that the wave functions be regular at the origin and damped at oo, we will construct 
the wave functions in the forbidden regions so that they fall off as r moves away from the 
turning points. Choosing 



Af 



N 



(32) 



The solutions in Region I require that the sign of the root ^o(^) — ^ko{r) be chosen so 
that the wave function vanishes at the origin. The radial wave functions in this region are 



^ , , NJm + S{r) - V{r) + E 
Gj{r) = —— 1 exp 

koir) 



dr'Cko{r') + hir')) 



(33) 
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and 



Frir) 



N 



X exp 



The solutions in Region III are 



^Jh{r)^m + Sir) - V{r) + E 

dr'{koir') + hir')) 



G 



iii[r) 



and 

Fniir) 



^^^^^N y/m + S{r)-V{r) + E 
ko{r) 



exp 



dr'Ckoir') + hir')) 



(34) 



(35) 



X exp 



\lh{r)^m + S{r) - V{r) + E 
dr'Cko{r') + hir')) 



-53(f;(r))-^o(r)-^i(r) 



(36) 



It should be pointed out that this is very similar to the WKB approximation to the Dirac 



equation described in 



53 



0|. In this previous work the WKB approximation is assumed 
to be in the form of a two-dimensional spinor of amplitude functions multiplied by the 
usual exponentiated phase function. Since the radial Dirac equations only determine two 
functions, it was necessary to make a choice for one of the amplitude functions that gives 
an upper-component wave function identical to that obtained here. The other two functions 
were obtained by solving two-dimensional matrix equations by introducing a dual set of 
spinors for a nonhermitian matrix and projecting to obtain scalar expressions for the phase 
function and the remaining amplitude function. However, this approach does not provide an 
expression for the first-order corrections to the phase that come from ^i(r) in the derivation 
given here. These contributions are necessary for a simple, smooth extrapolation of the wave 
functions over the regions where the WKB approximation does not converge. 



III. MODEL CALCULATIONS 

Our model consists of a constituent quark bound to an infinitely heavy di-quark and is 
represented by the Dirac hamiltonian for a particle in a scalar field S{r) and a time-like 
vector field V{r). 

H = a-p + l3{m + S{r)) + V{r). (37) 
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TABLE I: Comparison of eigenenergies calculated with RKF and WKB methods for a few selected 
states. 



11 




^RKF (GeV) ^wKB (GeV) ^rkf 


— -C/WKB 


1 


-1 


0.687666 


0.6895116 


■0.001846 


60 


5 


6.697930 


6.6979039 


0.000026 


196 


-9 


11.996600 


11.9965140 


0.000086 


20 


13 


4.353070 


4.3530336 


0.000036 


130 


-22 


10.051700 


10.0516770 


0.000023 



The scalar potential is a linear confining potential given by 

S{r) = br, 6 = 0.18GeV^ (38) 

In our model, the vector potential is provided by a vector color Coulomb potential. We will 
calculate for the case where the vector color Coulomb potential is absent, that is V{r) = 0, 
and where the vector potential is the simple static Coulomb potential 

V{r) = -^, (39) 

with /3 = 0.4. For convenience, the mass has been chosen to be m = 0. For these potentials it 
is easily shown that the Dirac WKB wave functions have the correct functional dependence 
near r = for all values of k. 

Note that we assume that the virtual photon only interacts with the light quark, and not 
with the infinitely heavy di-quark. This still allows us to gain qualitative insight into the 
issues of scaling and duality, but makes a direct comparison of numerical results from our 
model to experimental data impossible: our model is much closer to electron scattering from 
a B meson than a proton. We would like to point out that therefore, the values obtained e.g. 
for the momentum transfer at the onset of scaling should not be compared to experimental 
values obtained in inclusive electron scattering from the proton. However, the question if 
the bound-bound and bound-free transitions lead to the same scaling curve can be discussed 
within this model. 

One measure of the quality of the WKB approximation is a comparison of the eigenen- 
ergies calculated with the RKF method and the WKB approximation. Table |l] shows the 
eigenergies for a few selected states. Since the WKB approximation should be most accurate 
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FIG. 2: Reduced radial wave functions for k = —22 and n = 130 for a linear confining potential. 
Upper and lower components are shown for direct integration of the Dirac equation using the 
Runge-Kutta-Fehlberg (RKF) method and in the WKB approximation. 

at large energies and angular momenta, it is not surprising that the largest difference is for 
the ground state (lowest positive energy state), although even here the difference is less that 
2 MeV. For the other states the difference is less that 0.1 MeV. 

Figure 121 shows the upper and lower component wave functions for (3 = with k, = —22, 
n = 130 and E = 10.051 GeV (the last state listed in Table H}. Wave functions are shown 
for both direct integration of the Dirac equation using the Runge-Kutta-Fehlberg (RKF) 
method and for the WKB approximation. The range in r is chosen to cover the region of 
appreciable overlap between this state and the lowest positive energy state. The differences 
between the two sets of wave functions can not be distinguished in this figure. This indicates 
that the Dirac WKB approximation is very good at these energies. 

Figure El shows the longitudinal and transverse response functions for three-momentum 
transfers of 10, 20 and 30 GeV as a function of the y-scaling variable. Each plot contains 
a comparison of the model with WKB wave functions to the response calculated assuming 



13 



TABLE II: Extrapolation of the peak value -Rmax and position ymax of the peak of the response 
functions for g — > cxo. AR is the difference in peak height between the bound- free and bound-bound 
calculations. 







bound-bound 


bound-free 








2/max 


-Rmax 




Rmax 


AR 


linear 


Rl 


0.102 


0.683 ± 0.005 


0.105 


0.685 


0.002 




Rt 


-0.106 


1.358 ± 0.008 


-0.105 


1.370 


0.012 


linear -|- Coulomb 


Rl 


0.150 


0.606 ± 0.003 


0.132 


0.614 


0.008 




Rt 


-0.108 


1.137 ± 0.005 


-0.132 


1.229 


0.092 



that the quark is ejected from the bound state to a plane-wave continuum state, the bound- 
free transition. This bound-free result is the model equivalent of the parton model. In 
addition the plots for g = 10 GeV also show the model response calculated with RKF wave 
functions. It is clear from these that the RKF and WKB responses are virtually identical, 
as would be expected from the comparison of wave functions in Fig. |21 This figure shows 
that model response functions approach the corresponding bound-free response functions as 
the momentum transfer increases and that at g = 30 GeV they are almost identical. This 
indicates that the bound-bound transition within a model with a linear confining potential 
is dual to the bound-free case. 

Figure m shows similar response functions for the case where the vector Coulomb potential 
is included. Again, the RKF and WKB response functions are very similar at g = 10 GeV. 
As the momentum transfer increases, the longitudinal model response for the bound-bound 
transition appears to be approaching the bound-free result from above. For the transverse 
response, the bound-bound response is substantially smaller than the bound-free. Although 
the bound-bound response does slowly approach the bound-free with increasing momentum 
transfer, it is still quite far from the bound-free at g = 30 GeV. Also note that for the 
both longitudinal and transverse responses the shapes of the model result and the bound- 
free results are slightly different and this does not seem to be changing with momentum 
transfer. This seems to indicate that the bound-bound and bound-free responses are not 
dual once the Coulomb potential is introduced. An extrapolation of the peak positions and 
values of the response functions has been performed fitting the bound-bound calculations at 
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FIG. 3: Longitudinal and transverse response functions for the linear confining potential (/3 = 0) 
at q =10, 20 and 30 GeV. At each momentum transfer the bound-bound response function is 
compared to the bound-free response. The corresponding responses using the RKF integration are 
shown for comparison. 
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FIG. 4: Longitudinal and transverse response functions for the linear confining plus Coulomb 
potentials (/3 = 0.4) at q =10, 20 and 30 GeV. At each momentum transfer the bound-bound 
response function is compared to the bound-free (bf) response. The corresponding responses using 
the RKF integration are shown for comparison. 
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10, 15, 20, 25 and 30 GeV to the form 

B C 

i?max(g) =A+- + -. (40) 

q 

An estimate of the extrapolation error was obtained by varying the number of points in 
the fit from 3 to 5. Asymptotic values for the bound-free case are obtained directly from a 
direct calculation of the responses in the limit q oo. The results of the extrapolation are 
shown in Table IHl For the case of the linear potential the positions and peak values of the 
bound-bound response functions do not very significantly from the asymptotic bound-free 
response functions. However, for the linear-plus-coulomb potential, there is a consistent 
shift in the peak position and the difference for the height of the transverse response is 
substantial. A similar situation has been reported elsewhere for the case where the scalar 
a„d vecto. pote,.ia. a. — to be ident.ca. >.ea. confi.ng p„ten.a. QQ. 

In this paper, we have presented a new method to solve the Dirac equation for scalar 
and time-like vector potentials in the WKB approximation. We have applied this method 
to calculating the response functions of a light quark bound to an infinitely heavy di-quark 
in the bound-bound transition at very high momentum transfers. This type of calculation 
is relevant for modeling quark hadron duality. We compared these results to the bound-free 
transition, and found that the responses scale to the same limit for just a scalar potential. 
The vector potential introduces small differences in the two scaling functions. 
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APPENDIX A: THE DIRAC WKB WAVE FUNCTIONS NEAR THE TURNING 
POINTS 

Since the WKB expansion does not converge near the classical turning points, the WKB 
wave functions near these points are not a good approximation to the exact wave functions 
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and are, in fact, singular at the turning points. This problem can be eliminated by matching 
the WKB wave functions at some distance from the turning points to an exact or approx- 
imate solution to the wave equations near the turning points. The usual approach is to 
approximate the solution near the turning points by assuming that the potential is roughly 
linear over the region where the WKB wave functions fail. This gives matching conditions 
that lead to the phase in (j^^ and the corresponding wave functions are Airy functions Ai{z). 
The matching to the WKB wave functions is achieved by using the limiting properties of 
the Airy function 

lim n^z^Ai (z) -^e"^ (Al) 

and 

lim TT^z^Ai {—z) cos ( C — t ) ' (^^2) 

z—>oo \ 4 / 

where ( = |z5 can be identified as the WKB phase function. The factors of z^ cancel the 

_ 1 

singularity coming form the factor of ^ that appears in the definition of the WKB wave 
functions and gives a smooth finite result at the classical turning points. 

A similar procedure can be used for the Dirac WKB wave functions once a few additional 
complications are dealt with. The primary problem is that the local wave vector fci(r) is 
singular at the turning points. This is a result of the factor of kQ{r) that appears in (|19|). It 
is easy to show that near the turning points ko(r) ~ |r — r±|3. Therefore, ki(r) ~ \r — r±\^^. 
Since this singularity is integrable, ^^i{r) will be finite at the turning points but will have 
infinite slope. This phase also changes sign at these points. One result of this is to cause 
the total phase of the WKB wave functions to have an additional zero near the classical 
turning point. This interferes with the use of the usual Airy function approximation to the 
wave functions near the turning point for the upper component wave function. 

For the lower component wave function, there is the additional complication that it 
contains an explicit factor of ki{r) and is therefore more singular at the turning point than 
can be cancelled by the Airy function approximation, provided that a function that smoothly 
extrapolates to the sine function can be constructed. There is an additional problem arising 
from the factor of k'Q{r)/ko{r) that appears in (j2(J|) . This behaves like |r — r±|~i near the 
turning points and cannot be cancelled by the Airy function solutions. 
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Consider the solution in Region II given by (PHj). Define the phases 



C/(r) = 


/ dr'ki{r') 

Jr 


(A3) 


Cf-(r) = 


pr 

/ dr'ki{r') 


(A4) 


Crir) = 






/ dr'ki(r') 

J r 


(A5) 


crir) = 


/ dr%{r') 


(A6) 




Jr+ 





where i = 0, 1. Equation (j2Hl) can then be written as 
Gn{r) = N- 



^ym + S{r) - 


- V{r) + E 






\/m + S{r) - 


- V{r) + E 



-sin(Co^^-(r)-^)sin(Cr-(r))) (A7) 

Near the lower classical turning point, cos (Co (^) ~ f ) can be replaced by the Airy function 
using ()A2|) and sin (Co — f ) can be replaced using the identity 

lim -K^z'^'^Ai' {-z) sin f C - - ) • (A8) 

Therefore, near the lower classical turning point in Region II we replace the WKB wave 
function (j2Hl) by 



^ / X \/m + Sir) - V{r) + E i i , 
Gii{r) - -TfTs^ vr^^3(Co'^-(r)) 

X [Ai (-^(Co''-(r))) cos(Cr-(r)) - , ^ (-^(Co''-(r))) sin(Cr-(r)) 

V ^Jz{C^'-{r)) 

(A9) 

This can be continued across the lower turning point into Region I using ()A1|) and the 
identity 

- lim n^z'^Ai' {z) -e""^ (AlO) 

2— >00 2 
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12 3 4 

r(GeV-^) 

FIG. 5: The WKB wave function is shown as sohd hne. The Airy function, which smoothly 
interpolates between the WKB wave functions in regions I and II, is shown as dash-dotted line. 



to give 



Gi{r) 



N- 7=== n^z^iCoir)) 



X Az (^(Co'(r))) cosh(C^(r)) + ^-^ Az' {z{C^{r))) sinh(C^(r)) 



(All) 



These wave functions are now continuous and smooth at the lower classical turning point, 
and match (j^H|) and (jHSI) away from the turning points. The solutions near the upper turning 
point are obtained by the replacement C,i^~{r) Ci^'^i^) i^i (!A9|1 and for Region III by the 
replacement C/(^) ~^ (i^^i^) iii (!A11|) . 

The corresponding expressions for the lower components can be obtained by using the 
expressions above in 

G'{r) + f G(r) 



^^""^ ~ (m + S{r) - V{r) + E) ' ^^^'^^ 
An example for the procedure outlined above is shown in Fig. El The solid line shows the 
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WKB wave function given by ()33|) and ()28|). which clearly becomes singular at the turning 
point. The dash-dotted line is calculated using ()A11|) and ()A9|) . The interpolation of the 
wave function using the Airy function clearly matches the WKB wave function within one 
wave length on either side of the classical turning point. Equations (j33|) and (j28|) can then 
be used to interpolate the WKB wave function through the region of the classical turning 
point where the WKB approximation is invalid. 
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